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Abstract 

In this paper we propose and analyze a class of simple Nystrom discretizations 
Q of the hypersingular integral equation for the Helmholtz problem on domains of the 

Q plane with smooth parametrizable boundary. The method depends on a parameter 

(related to the staggering of two underlying grids) and we show that two choices 
of this parameter produce convergent methods of order two, while all other stable 
methods provide methods of order one. Convergence is shown for the density (in 
^ uniform norm) and for the potential postprocessing of the solution. Some numerical 

^ experiments are given to illustrate the performance of the method. 

•S 

3 1 Introduction 

In this paper we propose and analyze a discretization of the hypersingular integral equa- 
tion for the Helmholtz equation on a smooth parametrizable simple curve F C M^: 

oo 

O Here ' is the Hankel function of the first kind and order zero, k is the wave number, d^^ 

^ is the normal derivative, and g is data on F. The method is based on some simple ideas: 

(a) The operator is first written as a bilinear integrodifferential form acting on periodic 
functions 



9,^9,(y)//«(A;|--y|)0(y)dF(y) = ^ on F. (1) 



Jo Jo 

-k'^ [ [ H^^\k\yi{s) - x{t)\)n{s) ■ n{t)j/j{s)ip{t)dsdt, 
Jo Jo 
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X = {xi,X2) being a regular parametrization of F, and n = (xg, —x'l) the outward 
pointing normal vector. 

(b) The principal part of the bilinear form (the one with the derivatives) is formally 
approximated with a nonconforming Petrov-Galerkin scheme, using piecewise con- 
stant functions on two different uniform grids with the same mesh-size h: {{j — ^)h} 
and {(j — ^ + e)h}. Since the derivatives of piecewise constant functions are lin- 
ear combinations of Dirac delta distributions, that part of the bilinear form is just 
discretized with the matrix 

Wi+ij+i - Wij+i + Wij - Wi+ij, 

where W.j = H^'\k\hl - hj\), bf := x((2 - ^ + 6)h), and b,- := x((j - l)h). 

(c) The second part of the bilinear form (which has a weakly singular logarithmic singu- 
larity) is discretized with the same Petrov-Galerkin scheme, using midpoint quadra- 
ture to approximate the resulting integrals: 

eH^,'\k\mt-in,\)nt.n, 

j-{i+\+e)h j-(j+\)h 

^k^ / ifJ^^(A;|x(s) -x(t)|)n(s) ■n(t)dsdt, 

where mf := x((z + e)h), nf := hn{{i + e)h) and nij and rij are similarly defined. 

(d) The right-hand side is tested with piecewise constant functions and then midpoint 
quadrature is applied to all the resulting integrals. 

As can be seen from the above formulas, this method leads to a very simple discretiza- 
tion of ([T]), requiring no assembly process, no additional numerical integration and no 
complicated data structures to handle the geometric data. 

The use of a two-grid Nystrom method for periodic logarithmic integral equations 
goes back to the work of Jukka Saranen, Ian Sloan and their collaborators [HI HDl [13]. It 
was then discovered that the values e = ±1/6 provide super convergent methods (of order 
two) and that the values e = ±1/2 lead to unstable discretizations. The idea was further 
exploited in [5], showing that the methods can be used on the weakly singular equations 
that appear in the Helmholtz equation. The present paper shows how to transfer the same 
kind of ideas (and, up to a point, the same type of analysis) to the hypersingular integral 
equation for the Helmholtz equation. The case of the Laplace hypersingular equation is 
included in the present analysis. 

We will show that this discretization of the hypersingular integral operator is stable in 
the norm for the underlying space of piecewise constant functions as long as £ 7^ ±1/2 
(the value e = is excluded as a possibility from the very beginning, since it leads to 
evaluation of the kernel functions on the singularity). We will also show that e = ±1/6 
define methods of order two and that this order is actually attained in a strong norm. 
The error analysis will be based on Fourier techniques [H [9l [11] combined with an already 
quite extensive library of asymptotic expansions developed by two of the authors of this 
paper with some other collaborators O [H [5] [6] . 
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In a forthcoming paper p] we will show how to combine this discretization method 
for the hypersingular equation with the original method [21 [9l [TOl [13] for the single layer 
operator and with straightforward Nystrom discretization of the double layer operator 
and its adjoint. This results in a compatible and straightforward-to-code fully discrete 
Calderon Calculus for the two dimensional Helmholtz equation on a finite number of 
disjoint smooth closed curves. This discretization set has a strong flavor to low order 
Finite Differences. This might make it an attractive option to build simple code for 
scattering problems, when the simultaneous use of several boundary integral operators is 
required. 

The paper is structured as follows. In Section [2] we present the method for a class 
of periodic hypersingular equations that include ([T]) after parametrization. The method 
is then reinterpreted as a non-conforming Petrov-Galerkin discretization with numerical 
quadrature. In Section[3| we introduce the functional frame for the analysis of the method, 
based on the theory of periodic pseudodifferential operators on periodic Sobolev spaces. In 
Section |4] we present the stability result of this paper in the form of an infimum-supremum 
condition. In Sections [5] and |6j we respectively give the consistency and convergence error 
estimates for the method. Section [7] contains some numerical experiments, while we have 
gathered in Appendix |X] the more technical proof of Proposition 12 



2 The equation and the method 

We are going to present the method for a hypersingular integral equation associated to 
the two dimensional equation on a simple smooth curve. The extension to a finite set of 
non-intersecting smooth curves is straightforward. Since everything can be expressed with 
periodic integral equations, and the analysis will be carried out at that level of generality, 
we start this work by presenting the equation in that language. 

2.1 A class of integrodifferential operators and its discretization 

We consider two logarithmic kernel functions: 

Veis,t) := ^Ms,t)\ogism^7T{s-t)) + Keis,t), i E {1,2}, (2) 
where Ai, A2, Ki, K2 G C°°(]R^) are 1-periodic in both variables and 

Ai{s,s) =ay^O. (3) 
We consider the associated integral operators 

V,^:= [ V,{-,tMt)dt, (4) 

and the operator 

W:=-DViD + V2, Bip:=ip'. (5) 

Let us consider the space of smooth 1-periodic functions V := {ip E C°°(R) : ip{l + • ) = 
if}. An important hypothesis is injectivity: 

<peV, W<p = =^ <p = 0. (6) 
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As we will see later on, this injectivity condition is enough to prove invertibility of W in 
a wide range of Sobolev spaces. Given a smooth 1-periodic function we look for ip such 
that 

Wif = g. (7) 

The discretization method uses four sets of discretization points. Let be a positive 
integer, h := and 

Si = {i-l)h, ti:=ih, Si+s ■= {i + ^ - l)h, U+e ■= {i + s)h, i e Z. (8) 

(We will comment on e shortly.) The discretization method looks for 

N-l 

((^0, ■ • • , ^N-i) e such that ^^3^1 = h9{U+e) 1 = 0,..., N-l, (9) 

j=0 

where 

+ h^V2{ti+e,tj). (10) 

Substitution of £ by £ + 1 produces the same method. The option e G Z is not practicable, 
since it leads to evaluations of the logarithmic kernels in their diagonal singularity. The 
method fore = 1/2 (eG 1/2 + Z) will not fit in our analysis, that relies on stability 
properties of an e— dependent discretization of logarithmic operators that is unstable for 
e = 1/2. (We will show numerical evidence that the value e = 1/2 is valid though.) All 
other methods will provide convergent schemes, with two superconvergent cases. Namely, 
we will see that for smooth enough solutions, we can prove: 

max \ipj — (p{tj) \ = 0{h), if e ^ IZ 

3 

(this excludes the non-practicable and unstable cases), and that 

max|(^j - (^(tj)| = C'(/i^), lie e ±1 + 1. 

These results will be proved as Theorems [3] and |4] respectively. 



2.2 A non-conforming Petrov-Galerkin method 

We next give some intuition on how to come up with the method (|9|)-(|l0|)- We can formally 
rewrite ([T]) in variational form 

ij'{s){Vi^'){s)ds+ I i,{s){Y2^){s)ds= I i^{s)g{s)ds. (11) 
Jo Jo 

Consider now the function Xi that arises from 1-periodization of the characteristic function 
of the interval (sj, Sj+i) = (tj — h/2,ti + h/2), that is, 

X.(l +■)=X^. = 1, = 0. (12) 



We similarly define the functions Xi+e by periodizing the characteristic functions of the 
intervals (sj+e, Sj+i+j) = (tj+e — /i/2,tj+e + h/2). The weak derivatives of these functions 
can be expressed through the use of Dirac delta distributions. In addition to under- 



standing Dirac deltas as periodic distributions (see Section 3.2), we will admit the action 
of Dirac deltas on any function that is continuous around the point where the delta is 
concentrated. At the present stage, we only need to consider the functionals 



{Ss^, if} := ip{si), {Ss^+^, ^} := ^isi+,), 



(13) 



acting on any 1-periodic function ip that is smooth in a neighborhood of Sj and Si+e- 
Admitting formally that Xi = ~ ^si+i cind x'i+e = 



we consider a non- 



conforming Petrov-Galerkin discretization of (11): 

N-1 



N 



where 



ifjXj such that ^ W^fipj 



j=0 



j=0 



g{s)ds, 



N 



ti+e 



(14) 



Note that 



^l('5i+l+e, ■Sj+l) — ^l('Si+£, Sj+i) — Vi(Sj+i+e, Sj) + V^i(Si+e, Sj) 



+ 



J. _L-ik ni J-H. 



I h_ I • h 



V2{s,t)dsdt. 



N-1 

E 

j=0 



V2{s,t)dsdt](pj 



(y2Vh){s)ds 



and that the leading term in W^j° can be understood as the action 



(15) 



(16) 



The method ([9])-(10) is recovered if we use midpoint integration for all integrals in (14) 
(p|. 



2.3 Relation to the Helmholtz equation 



Let X : 



— > 



be a smooth 1-periodic function such that |x'(s)| ^ for all s, and 



x(s) 7^ x(t) if |s — tl < 1. The range of x is then a smooth closed curve T in the plane. 
Let n{t) be the outward pointing normal vector at x(t) with |n(t)| = |x'()f:)|. Given a 
periodic function ip, we define 



\r 9 z ^ f/(z) 



(1)^ 



|z-y|) 



n(t) Lp{t)dt 



4 -/o - " ^ ' " " y=x(t) 

"ik f^^ju,,, , , ,, (z — x(t)) ■ n(t) 

iff )(A;|z - x(t)l) |z-x(t)| ^^^^"^^^ 



(17) 
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where i^Q^^ and h[^^ are the Hankel functions of the first kind and orders and 1 respec- 
tively. The function U is an outgoing solution of the Helmholtz equation 

dU 

AU + k'^U = inM^r, _ ^kU = oir-^^^) as r = Izl ^ cx), (18) 

or 

with the asymptotic limit (the Sommerfeld radiation condition) holding uniformly in all 
directions: U is the double layer potential with (parametrized) density ip (see |7l Section 
2.1]). The double layer potential is discontinuous across T but its normal derivative on T 
coincides from both sides. If we define 

iWip)is) :=-(Vf/(x(s))-n(.), 

then, 

{Wv){s) = -^A^'i/«(fc|x(.)-x(t)|)^'(t)dt 



H!,'\k\^{s) - x(t)|)n(s) ■ n(t) ^(t)dl (19) 

Jo 



This is just the parametrized form of a well known formula: see [12], Theorem 3.3.22] or P, 
Exercise 9.6]. Using the asymptotic behavior of Hankel functions close to the singularity, 
it can be shown that the weakly singular kernel ifo^'*(A;|x(s) — x()f:)|) can be decomposed 
as in ([2]) with Ai(s, s) = i/2 (see ([3])). Therefore, the operator W in (19) fits in the frame 
of (§. 



The operator W satisfies the injectivity condition ^ if and only ifk^ is not a Neumann 
eigenvalue of the Laplace operator in the interior of T [7, Section 2.1]. In those cases. 



the solution of the interior-exterior Helmholtz equation (18) with Neumann boundary 



condition d,yU = f can be represented with the double layer ansatz (17), (p being a 
solution of equation ([7]) with g := — / o x. If we then apply the discretization 
(10), we can construct a fully discrete potential representation, when we substitute (f by 
iph '■= J2j^=o fjXj ill ( fl7| ) and apply midpoint integration: 



3 Functional frame 



3.1 Asymptotics of hypersingular operators 

Consider the space V of periodic C°° complex valued functions of one variable, endowed 
with the metric that imposes uniform convergence of all derivatives [11, Section 5.2]. A 
periodic distribution is an element of V, the dual space of V. Given u G V, we consider 
its Fourier coefficients 

u{m) := {u,(j)m)v'xv, 4>m{t) ■■= exp(27rzmt), m G Z. (21) 
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The periodic Sobolev space of order r G M is 



■.= {ueV' : \\u\\r<oo}, where \\u\\l := \u{0)\^ + J^H^^'l'^Ml^ ■ (22) 

(From here on, the symbol X^m^o I'sfers to a sum over all integers except zero.) An 
extensive treatment of these spaces can be found in the monograph [TT]. Let us just 
mention that C if for p > q, with dense and compact injection. Also, can 
be identified with the space of 1-periodic functions that are locally square integrable or 
equivalently, with the 1-periodization of L'^{0, 1). 

We say that an operator L : V ^ V is a periodic pseudodifferential operator of order 
n, and we write for short L G S^n), when L : — )■ if"" for all r. It then follows 
from [m Paragraph 7.6.1] that the logarithmic operators Q can be extended to act on 
all periodic distributions and, consequently, so can W. Moreover, Vi,V2 G £{—!) and 
W G ^(1). 

A first group of pseudodifferential operators that we will use extensively is that of 
multiplication operators. Given a G we define the operator a G £{0) by au := au. The 
periodic Hilbert transform 

Hm := sign(m)u(m)0m (23) 

is clearly a periodic pseudodifferential operator of order zero. We also consider the oper- 
ators for n G Z: 

D„M := ^(27r«m)"M(m)0^, D„ G £(n). (24) 

m^O 

It is easy to note that Di = D is the differentiation operator and D_i is a weak form of 
the following antidifferentiation operator 

(D_im)(s) = [ {u{t) - u{0))dt \Ju G 
Jo 

In the next lemma we collect some elementary properties of these operators. 



Lemma 1. The following properties of the operators D„ in (24) hold: 

(a) DqU = u — n(0) for all u. 

(b) For all n and u, D„„D„u = u — {i(0). 

(c) For all n,m, D^D^ = Bn+m- 

(d) For all n, D„H = HD„. 

(e) For all a E V, Dia = aDi + a' . 

The following results show how logarithmic operators and the hypersingular operators 
W can be represented up to operators of arbitrarily negative order as a linear combination 
of compositions of the simple operators given above. 
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Proposition 1. Let 



Vu:=—I A{- ,t)\og{sm^{n{- -t)))u{t)dt+ I K{- ,t) u{t)dt, (25) 
27r Jo Jo 

where A,K& C°°{M?) are 1-periodic in each variable. Then there exists a sequence 
{(^n}n>i C V such that for all M, 

M-1 

V = ^ a„HD_„ + Km, Km G £{-M). (26) 

n=l 

Moreover ai(s) = A{s, s). 

Proof. See |5l Proposition A.l] or similar arguments in pT| Chapter 6]. □ 

Proposition 2. LetW be the operator in ([s]). Then, there exists a sequence {bn}n>o C T) 
such that for all M , 

M-1 

W = -«HDi + bn'i^^-n + Km, Km G £{-M), 

n=0 

where 7^ a G C zs the constant in ([s]). 

Proof. Tliis is just a direct consequence of Proposition [T] and Lemma [TJ □ 

Proposition 3. Hypothesis ^ implies invertibility ofW : H'^ — > H^^^ for all r. 

Proof. Consider tlie lowest order expansion of Proposition |2| namely W = — aHD + Kq 
with Kq G £{0). It is clear that — aHD : — )■ H^^^ is Fredholm of index zero, and 
therefore so is W. 

Let now u G be such that Wu = 0. Applying that ¥Pv = f — v{0) for all v, and 
Lemma [l](b) and (d), it follows that 

= HD_iWm = -au + w, where w := au{0) + HD_iKom G H''^'^. 

This means that u G for all s and therefore u E V. Hypothesis ^ implies then that 
u = 0. Therefore W : — )■ H^~^ is injective and, by the Fredholm Alternative, it is 
invertible. □ 



1 /.I 



3.2 Variational formulation of the discrete method 

For a fixed t G M we can define the Dirac delta distribution 6t by its action on elements 
of V, {St,4>)v'xv = <P{t)- Using the Sobolev embedding theorem [HI Lemma 5.3.2], we 
can prove that St G for all r < —1/2. However, this does not allow us to apply 
the Dirac delta to functions that are piecewise smooth on points where they do not have 
jumps. If M is a 1-periodic function that is continuous in a neighborhood of t, we will write 
{6t,u} = u{t). Note that in general this is not a duality product V x V or if x H^^ . 
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With this definition, we can admit the Dirac deltas (5^. and 5s in (13), as well as formula 

Let Pq be the space of constant functions. We then introduce three iV-dimensional 
spaces: 

Sh := span{xi : i = 0, . . . , N - 1} = {uh e H° : Uh\(s,,s,+,) G Pq Wi}, 

.,iV-l} = {iz,, G//° : uj,\(s 



s, 



h,e 
r,-l 



span{x 



i+£ 



0, 



span{5s, : ^ = 0, . . . , - 1} = {m'^^ : Uh E Sh} © span {dh}, 



where 



Af-l 



(27) 



(28) 



N-l 

Finally, we consider the discrete operators 

N-l 

j=0 j=0 

that are well defined for all periodic functions that are continuous around U{tj} = hZ 
and U{tj+e} = h{€ + Z) respectively. In particular, we can apply Qj^^ to elements of Sh 
and Qf^l to elements of Sh,e- 

Proposition 4. Let £ ^ Z and let g be continuous in a neighborhood of h{e + Z). Then 
the discrete variational problem 



IS 



find iph G Sh such that 

{^;, Vi^'J + {Q;:>, V2Q^ V4 = {Qhli'h, g} ^i^k e Sh,e, 

equivalent to looking for iph = Ylj'=o fjXjj where equations (j9])-(10) are satisfied. 



(29) 



Proof. Given that both Sh and Sh,e are A^- dimensional, the problem (29) can be reduced 
to a X linear system, after choosing a basis for each of the spaces. The result follows 
then from several simple observations. First of all x'j = — ^s^+i and Xi+e = 5s,+, — ^Si+i+^ 
and {5s,+,,Vi5sJ = Vi{si+e,Sj). Also 

Finally {Qh\Xi+e^g} = hg {ti+^ ) . □ 



The non- conforming Petrov-Galerkin discretization of (11) given in (14)-(15) is equiv- 
alent to the discrete variational problem 



find iph G Sh such that 

{ijl Vi/J + iiJh, y2^h) = ii^h, g) ^i^h e Sh, 



(30) 



where (m, v) = u{t)v{t)dt. In the sequel ( ■ , ■ ) will be used to denote this bilinear form 
m (so that ( u.,u) — ll'^^llo) and its extension to a duality product H ^ x H^, so that 
for any r G M. 

= sup i;Gif. (31) 



\u\ 



We will always take adjoints with respect to this bilinear form, thus avoiding conjugation. 
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4 Stability analysis via an inf-sup condition 



Consider now the bilinear form Wh '■ Sh^s x 5*/! — ?■ C associated to the problem (29 ), namely 

Whi^ph, Vh) ■■= {^/^l + {Qhi^h, VsQ^ V4- (32) 
The aim of this section is the proof of the following result, that in particular implies that 



problem (29) (and by Proposition |4] also the method (|9])-(10)) has a unique solution for 
small enough h. 

Theorem 1 (Stability). There exist and /^^ > such that for h < Hq, 

\Wh{%l)h,(ph)\ 



sup 



> /3e\\^h\\o Vv9/j E Sh- 



4.1 Stability of the non-conforming PG method 



We start by considering the bilinear form associated to problem (30) 



(33) 



the operator Aip := WD_iy9 + ^(0)W1 and its adjoint A*. Note that if iph ^ Sh,e, then 
A*iph is a smooth function except at the discontinuity points of iph- 

Lemma 2. 

wl{iJh,<Ph) = Wh,-^*iph} + i^hio){i,A*iph) '^iph e Sh,s, e Sh- 

Proof. A direct computation, using the fact that Dl^^ = — D_i and Lemma [l|^b) shows 
that 

A*^' = -D_iW*^ + W*^(0) = V*V^' - \^'(0) - D_iV;^ + V^(0). (34) 
Noticing that {y^'^, 1} = (y?'^, 1) = 0, it follows from ( [34) ) and Lemma [TJ^b) that 

M,A*v^,} = {y.'^vt^a + lD-iv'Lv;^,) 

= ivi^i'^'h} + {vh.yi'^h) - Vh{m,yi'^h) = <{'4^h,vh) - ^ft,(o)v*^;(o). 



At the same time, integrating in (34), it follows that (1,A*'?/') = V2'?/'(0), which finishes 
the proof. □ 

Lemma 3. \\dh — < vr/i. 

Proof. Since for all u G H^, 



'I - dh,u) 



„i N-l 
/ M(t)dt - /i V m(s 



N-l 
j=0 



h 

J 2 



u{t)(M — hu{sj 



< 



h 



N-l 



h 



u{t)\dt < — ||m'||o ^ 7r/i||n||i, 



j=0 "'=3-2 



the result follows from (31) 



□ 
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Lemma 4. For e ^ \L, there exist positive constants (3^ and Ce such that 



sup 



|{5h,A*V'h}| 



and 



1(1, A*v^O - {4, A*^,}| < CMli^hh v^, G 



(35) 



(36) 



Proof. Using Proposition |2] (asymptotic expansion of W), it is simple to see that A = 
— aH + K, where K G 8{—l). This places us in the hypotheses of [2i Proposition 8], which 
proves (35). For the second estimate, we write A* = — aH + K* and decompose 

7V-1 

(1, A>,,) - {4, A>4 = -a((l,H^/,) - h Y,{^^h){sj)) + (1 - dh^KiPh). 

3=0 

The first term can be bounded using [21 Lemma 3] 



N-l 
j=0 



yiph G She- 



The second one follows from Lemma [3| 

|(l-4,KV^h)| < ||l-4||-i||K*^/.||i <vr/i||K*||^,o^^,i||^;,||o. 

Lemma 5. There exist two positive constants Ci,C2 such that 

CiWh + '?h{0) 411-1 < Ibftllo < C2\\^'h + ^h(O) 411-1 V(/3;, G Sh- 
Proof. The first bound is a simple consequence of the following inequalities 

llv^hll-i < 27r||(y9ft||o, \^h{0)\ < WVhWo, II4II-1 < 1 + tt/i 



□ 



(the last inequality follows from Lemma |3|) Note now that the operator Sh S^^ ^ given 
by (ph ^ ip'h + '?h{0)dh is injective. Its inverse is S^'^ 3 6h ^ D-i^/^ + Sh{0){l - D_i4). 
Then the second inequality of the statement follows from the fact that 



h 
2 ' 



||D_A||o < ^\M-u Mo)\ < \M-u 111 - D_i4||o < 1 + ^11411-1 < 1 + 

where we have applied Lemma |3] again. □ 
Proposition 5. There exist ho and Pe > such that for h < Hq, 

\w°hi'^h,<fh)\ 



sup 



> /3e||v5/i||o Vv^ft G Sh- 
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Proof. By Lemma [2| we can write 

wliiPh, ^h) = Wh + MO)dh, A*^4 + ^,,(0) ((1, A*V^fc) - {4, A*V4 
Therefore, by Lemma |4| it follows that 

sup TT—r. > sup — — Csh\i^h[0)\ 

0^Vh6Sh,s Wh\\0 Oy^i'heSh,e WWhWo 

> /3e\\(p'h + (ph{0)dh\\-l - Ceh\\(ph\\o 

> {l^sC2^ - Ceh) W^hh e Sh, 

where we have applied Lemma |5] in the last inequality. □ 

4.2 A perturbation argument 

Consider now the quadrature error 

E!j-= I V2{s,t)AsAt-h^V2{ti+e,tj), 1,3 eZ, (37) 

where := Sj+i+e) x (sj, s^+i). Let be the N x N matrix whose entries are the 
values Efj for i, j = 0, . . . , N — 1. In the sequel \E\p will denote the p-norm of the matrix 
E (for p e {1, 2, 00}) and |-E|Frob will denote its Frobenius norm. 

Proposition 6. 

K^h, V2v?h) - {Qh^^i'h,y2Qh^fh}\ < h''^\E^\2\\^ph\\o\\fh\\o ^i^h e Sh,e<^h e St. 

Proof. If we decompose iph = S^lo^ i'iXi+e and iph = J2f=o VjXji it is easy to see that 

7V-1 
i,j=0 

The result is then straightforward noticing that ||'?/'h||o = ^^^^K'^o, • • • , V'Af-i)! and ||v5/i||o = 
h^^'^K^Po, . . . , v97v-i)|, where | • | is the Euclidean norm in C^. □ 

In order to simplify some forthcoming arguments, let us restrict (without loss of gen- 
erality) e to be in [—1/2,1/2] \ {0} (the restriction e 7^ ±1/2 is not needed for these 
arguments). 

Lemma 6. There exists Ce such that for all h 

\E^<Ceh^\\ogh\ Vz,j. 
Moreover diverges like log \e\~^ as \e\ — )■ 0. 
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Proof. Since Efj_^j^ = E^j, we can choose |« — j| < N/2 and then 



\s-t\<h + h{^ + \e\) < |/i + I < |, {s, t) e (38) 

as long as h < 1/6. (This is not a restriction, since we are only missing the values 
1 < < 5 that can be incorporated by modifying the constants in the final bound.) Also 



\£\h< \ti+s-tj\ < |. 



(39) 



Because of the form of the kernel function V2 (see (|2])), in the diagonal strip D := {{s, t) : 
|s — t| < 3/4}, we can bound 



\V2{s,t)\<C,\og\s-t\-' + C,, {s,t)eD. 



(40) 



Therefore, by (39), 



\V2iU+,,tj)\<C,\og\E\-' + C^\ogh-' + C2, Vz,j 



h/2, then by (40) 



[ V2{s,t)dsdt < Cih^ {\og2 + log h-^) + C2h^. 



(41) 



The choice of indices \i — j\ < N/2 ensures that Q^- C D. If dist^Qlj, {(s,s) : s G M}) > 



(42) 



If, on the other hand, dist(Qfj-, {(s, s) : s G M}) < h/2, a simple geometric argument 
shows that Qfj C {{s,t) : t G (sj,Sj+i), \s — t\ < ch}, where c := v^+1/2 < e. Therefore 

r pch 

/ V2{s,t)dsdt <2Cih \ogu-^du + C2h^ = 2Cich'^{l-\ogc + \ogh-^) + C2h^. (43) 
JqI, Jo 



We can now gather the bounds (41 ), (42 ) and (43 ), rearrange terms and take upper bounds 
to prove the result. 



□ 



Lemma 7. There exists C independent of e such that for all h 



\EU<C 



- for i,j such that 2 < \i — j\ < 



N 
~2' 



Proof. If (s,t) G Qfj, 2 < |z - j| < A^/2, and \e\ < 1/2, then 

\s-t\> \ti+e-tj\-\t-tj\-\s-t,+e\ >h\i-j\-{\e\ + l)h> \h\i-j\. (44) 

Recall first the definition of D := {(s,t) : |s — t| < 3/4} given in the proof of Lemma |6} 
Taking derivatives of the kernel function V2, we can write 



''s,t) = ^^ + B2{s,t), 



\s-t\' 



dt^ 



S,t) = ^^_^2 +B4{s,t), 
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where for I G {1,2,3,4} the functions 5^ G C^{D) are bounded. Using the error bound 
for the midpoint formula in two variables, we can estimate 



24 (.,t)eQ|: V ds^ 
( 

24 V {s,i)eQf, 



|i — V 4 ^ K ~ Jp 



where we have applied (44) and the upper bound \i — j\h < 1/2. 
Lemma 8. There exists such that for all h, |-E^|2 < C^h^^'^. 



□ 



Proof. We first decompose the matrix E'^ = E^^^^ + E^^, where E^^^^ gathers all tridiagonal 
terms (modulo A^) of E'^ 



Etj, < 1 or 

0, otherwise. 



Using Lemma |6] and the fact that -Efrid only three non- vanishing elements in each row 
and column, it is easy to estimate |-Etridli + l-^tridloo < 3Ceh'^\logh\. Therefore, by the 
Riesz-Thorin theorem 



EU2 < lEUY'lEU'l' < 3C,/i2|log/i| 



(45) 



On the other hand, we can estimate the off-diagonal terms using Lemma [7| (recall that 
we can move indices so that |« — j| < N/2) 



N-l 



\e: 



■e 12 
off 1 2 



< \E: 



e |2 
off iFrob 



< 



i=0 2<\i-j\<N/2 



\E: 



•e |2 
ij\ 



N-l 



< 



1 



o|4 



00 ^ 



i=0 2<|j-j|<oo 



k=2 



Gathering (45) and (46), the result follows. 



(46) 



□ 



Proof of Theorem [1} Note that 

As a direct consequence of Proposition [6] and Lemma |8| we can bound 

\{^h,y2V>h) - {QH!e^h,y2Qh'v>h}\ < C^h'^'UhhyhWo e Sh,e, v>h e Sf,. (47) 

The proof is thus a simple consequence of this bound and Proposition [5] □ 
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5 Consistency error analysis 

The analysis of the consistency error is based in the careful use of estimates for quadrature 
error and the combination of asymptotic expansions of discrete and continuous operators. 
We start this section with some technical results that will be needed in the sequel. 

5.1 Estimates for quadrature error 

Lemma 9. The following bounds hold for all h and all e: 

(a) \{Ql^^i)h,u) - {il)h,u)\ < \{'KhY\\tljh\\o\\u\\2 for all iph e Sh,e and u e H'^. 

(b) \{Q^^^'iph,u)\ < llV^/illodl^^llo + 7r/i||M||i) for all i^h e Sh,e and u e . 

Proof. Using Taylor expansions, it is easy to prove the following well-known bound for 
the midpoint formula 



u{t)dt — hu{c) 



< 



\u{t)\dt, 



from where 

\{Qhl'^h,u) - {iph,u)\ 



N-l 



< 



\Mt)\\u"{t)\dt< 



•h\\o{2-K) \\u\\2. 







(4J 



□ 



This proves (a). To prove (b) we proceed similarly, showing first that 

\{Qlli)h,u) - {iph,u)\ < 7r/i||^/,||o||M||i \/i)h E Sh,e, u e H^, 
and then applying the inverse triangle inequality. 
Lemma 10. There exists Ce such that 

\{Q-H^e^h.y2Q-H^u]\ < C,||^^||o(||n||o + h\\u\\^) V^^ G Sh,e, U E H\ 

Proof. Let Uh := J2^=o ^ Sh be the midpoint interpolate of u onto Sh- A direct 



estimate shows that 



\uh\\o < \\u\\o + \\u - M/,||o < ||m||o + ;;7f ||M||i- 



(49) 



On the other hand, since Q/^^m = Qf^^Uh, it follows from (47) that 



< (llValli^o^^i^o + C,h^/^)\\^h\\o\\uh\\o- 



Applying (49), the result follows. 



□ 
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5.2 Discrete operators and expansions 

The truncation operator for the Fourier series 

FhU := ^ u{m)(pm where Aat := {m : -N/2 <m< N/2} 

meA]v 

gives optimal approximation properties in all Sobolev norms [11, Theorem 8.2.1] 

\\FhU-u\\s < {V2hy-'\\u\\r. r>s. (50) 

We can also define a discretization operator onto Sh based on matching the central Fourier 
coefficients 

DhU G Sh such that Dhu{m) = u{m) Vm G Aat. 

This operator is based on a class of spline-trigonometric projectors introduced in [T]. Here 
we will use it as introduced in |6j. The following property 

Ql\,2Fh^ = '0Dh (51) 

is a consequence of [2| Lemma 5]. 

Consider the 1-periodic functions such that {—lYilB^ restricted to (0, 1) is equal 
to the Bernoulli polynomial of degree i for all i. Consider also Cj := HB^. By comparing 
their Fourier coefficients [21 Section 3], it is easy to prove that 

Ci(t) = -2^1og(4sin2(7rt)) and therefore C^{±1) = 0. (52) 

Note that ±1/6 + Z are the only zeros of (7^. 

Proposition 7. Let ai,a2 G V and 

V := aiHD_i + azHD^a + Kg where K3 G ^(-3). 

Let then Li := ai and L2 := aiD — a2, and consider the operators 

RhU := Yu-YQl'FhU + hC^{- /h)UFhU, 

ThU := Yu-YQl\^^FhU + hC^{- /h + l/2)UFhU + h''C^{- /h + l/2)L2FhU. 
Then 

\\Rhu\\o + h\\Rhu\\i < Ch^\\u\\i W e H\ (53) 
||T;,m||o + /^||T/,^/||i < Ch^Wuh ^ueH\ (54) 

Proof. It is a direct consequence of [2|, Proposition 16]. □ 

Proposition 8. Let EhU := u - Dhu + hB^{- /h + l/2)Fhu'. Then 

WEhuW^ + h\\Ehu\\^ < Ch^uh Vm G 

Proof. The bound for Hi^ft-uHo is given in |5l Proposition 1]. The bound can be obtained 
with similar arguments (see the proof of [21 Proposition 16]). □ 
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5.3 Consistency error 



In the definition of the bihnear form (32), we only admitted discrete arguments. In this 



section we will admit a continuous second argument. The definition is equally valid. 



Proposition 9. Let (ph be the solution of (29) with right-hand side g = W(p. Then 

\wh{iph,^h - V^) + hC^{e){ilJh,aip)\ < Ch"^ WiphWoWv^h ^^/i ^ Sh,e, 
where a{s) := ^2(5,5). 
Proof. Note first that by definition of (fh 

w,{^,, - ^) = {Qj;^^^,, W^) - {^l Vi^') - {Q^_^^,, V2Q;: V} (55) 
= -iQ^l^h, DViD^) + (^,, DViDv;) 



■■■■Ti 



+ {Q-,l^k, y2Ql\Fuv - v)} + {Qt^h. V2(^ - Q-^'Fnv)} . 



■T2 



■T3 



In order to estimate Ti, we apply Lemma [9[a) with u = DViDip and note that DViD G 
S{1), to obtain 

^/illollv^lls- (56) 



To estimate T2, we apply Lemma 10 and the approximation properties of (50), so that 
IT2I < UhUWFhV - '^llo + 7rh\\Fhy^ - <^||i) < Ch^^h\\oM2. (57) 

To bound T3 we will apply Proposition [T] to the operator V2 (see Proposition [T]). It is 
simple to verify that 

{QHi^h,C,i-/h)u} = C,is)iQj;^^^,,u) VV^, G S,,,. 

Then, by Proposition [7| and denoting a(s) := ^2(3, s), 

= -hC_^{e){i)h,aip) 

+ hC^{e){iJh - a^) + hC^{e){Q^^h. - Fh^)) + {Qlli^h, Rh^) . 

V V V 

=:T3i =:r32 =;T'33 

Applying Lemma [9]^a) with m = a we can easily bound 

\Tzi\ < CeWMUvU (58) 
while Lemma [9|^b) applied to m = a((/? — -FXv') yields 

IT32I < h\\^ljh\\o{\\a{^ - Fhv)\\Q + h\\a{^ - Fhv)\\{) < Ch''\\ilJh\\o\M2. (59) 



Finally, we apply Lemma |9j(b) again, using the bound for R^ip provided by (53), which 
yields 



IT33I < Uh\\o{\\Rhiph + T^h\\Rh'p\\i) < C/i"||^/.||o||^||i. 



(60) 
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Inequalities (58)-(60) imply that 

IT3 + hC,{e){iJh, a^)\ < Chm,uUM2. 



Carrying (56), (57), and (61) to (|55|) the result follows. 



(61) 

□ 



Proposition 10. Let a he the constant in (jsj). Then 



Proof. Using (51), we can write 



wnii^h. ^ - D,^) = Vi(^' - Qj;\/,F,^')} + {Q^^^^Ph, V2Q, ^(^ - D,^)} . (62) 



=:T4 





=:T5 



To estimate T4 we will use Proposition [T] applied to Vi (see Proposition [T]). An easy 
computation shows that 

C,{ ■ /h + 1/2) u} = C,ie) iijl u) ^ij, e SH,e. (63) 

Let Li and L2 be the differential operators associated to the expansion of Vi in Proposition 



7| (note that Liu = au). By ( p33| ) and Proposition |7| it follows that 
T4 ^ 



+ hC,iE)ifMif' - F,ip')) -h'C,ie)i^PlUF,ify+i^lT^') . (64) 

:T43 



^/ 

--■■Til 



'■■T42 



Using (50) we can easily bound 



1^411 = \hC,{e){^PH, DLi(v;' - Fh^')\ < C/i||V^^||o||</^' - Fh^'\\, < C'hm^hhMU (65) 



and 



IT42I = \h^C,ie)i^,,,I)UF,,if')\ < Ch'W^Mvl 



Similarly, using the bound for given by (54), we estimate 

1^43! < ||V^/i||o||D7),V9'||o < C/l^||V^/i||o||(/3||3- 



Taking (65)- (67) to (64) we have proved that 

IT4 + haC,{e){ijiy,')\ < Ch^MoM3. 



We next estimate the term T5 in (62). Note that 



Qh\M . /h + l/2)u) = B,{1/2)Q~,\ = 0, 



(66) 
(67) 

(68) 

(69) 



by the coincidence of B^i with the Bernoulli polynomial of first degree in (0, 1). Therefore, 
using Proposition [s] and Lemma [lO||^b) it follows that 

l^sl = \{Qf,l^h,y2Qh'EHip}\ < Cs\\Mo{\\Ehip\\o + h\\Ehip\\i) < C',h^\\iJhWoM2. (70) 



The collection of (62), (68) and (70) proves the result 



□ 
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Corollary 1. Let (fh be the solution of (29) with right-hand side g = W(f. Then 

\whii'h,'^h - Dhip) + hC_^{e){i)h,aip - a(^")\ < Ch'^WiphWoW^h '^i'h e Sh,e, (71) 
and 

\wh{^h,^h- Dh^p) + hC^{e){Q^^^i)h,aLp - aip")\ < Ch'^\\ilJh\\o\\^h ^i^h & Sh,e, (72) 
where a(s) = ^2(5, s) and a is the constant in 



Proof. The first bound is a straightforward consequence of Propositions [9| and [lOj The 
bound (72) can be derived from the first and (48), although it has already been implicitly 
given in the proofs above. □ 

Proposition 11 (Zero order asymptotics). The following reduced estimate holds: 

\whii'h,Vh- DhV)\ < ChWiJhWoW^h yi'h^Sh,e- (73) 

Proof. If we go back to the notation of the proofs of Propositions [9] and 10, it is clear 
from (|57l), (|6l|, and (|70|) that 

IT2I + IT3I + IT5I < C/i||^^||o||<^||2 e Sj;}. 

Using (48) instead of Lemma |9](a), it is also simple to bound 

|Ti| <C/i||^,||o||</^||2 VV^, G5,j. 

For the operator T^m := 'Vu — 'VQ^\^^FhU + hC_i{ ■ /h+l/2)LiFhU, we can bound ||T^n||o + 
/i||T^m||i < C/i^||m||i [2, Proposition 16]. Using this bound instead of Proposition [7| we 
can prove that 

IT4I < C/i||^/.;,||o||<^||2 \fA e S^^l 
This finishes the proof. □ 

In order to set up clearly the precise formulas of the second term in the asymptotic 
expansion of Whi'iph, fh — Dhip), we need to consider the first two terms in the expansions 
of Vi and V2 given by Proposition [ij 

V2 = aHD_i + bHD_2 + K3, K3G^(-3), 
Vi = aHD_i + cHD_2 + J3, J3 G ^(-3). 



Proposition 12 (Second order asymptotics). Let (fh be the solution of (29) with right- 
hand side g = Wip. Let 



pe 



C^{e) (piaD - c)D - (aD - b)) + ^(D^ViD + VaD^ 



Then 



1=1 

Proof. See Appendix |X] 



□ 
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6 Convergence theorems 



The results of Sections [4] and [s] give a first simple estimate of the convergence of the 
method, showing that when e = ±1/6, the solution superconverges to the projection Dhf. 
We first recall that [5], Formula (5)] 



\Dh<f - LfWs < CK" 



m 



s < r < 1, 



s < 1/2. 



(74) 



Theorem 2. Let iph be the solution of (29) with right-hand side g = Wip and e ^ 

(75) 



Then 
Moreover, 



WVh - T^hVh < Ch^iph ife e {-1/6, 1/6}. 
Proof. Using Theorem [T| and (73), we can prove that 

\whiiph,'fh - Dh'f)\ 



(76) 



l^eW^h - Dh^Wo < sup 



<chyh 



Applying (74), this proves (75). The superconvergence bound (76) follows from Corollary 
[l| (note that ^^(±1/6) = 0) and Theorem |l| □ 

The superconvergence estimate can be first exploited with a postprocessing of the 
solution: given v smooth enough we approximate 

N-l 



/ (p{t)v{t)dt ^ hy^'Pj^i^j) = iQh^Vh,v). 
Jo ,=0 



This includes the fully discrete double layer potential (20) to approximate (17). 



Corollary 2. Let (ph he the solution of (29) with right-hand side g = W^p, and e e 
{-1/6,1/6}. Then 

Proof. Using Lemma |9Va) (with e = 0), we can easily bound 



< Cih^\\>ph\\Q\\v\\2 + W^h - -D/iV^lloll^^llo + \\Dh^ - p:>\\-i\\v\\i 

< C2h^ (\\'P>\\2\\v\\2 + llv^llsll^llo + IIV'lllll^lll), 



by (17§, (17§, and 

We next introduce the interpolation operator 

N-l 

hu := ^u{tj)xj. 

3=0 



□ 
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The Sobolev embedding theorem \TD; Lemma 5.3.2] and Proposition |8] show that 

max|M(t,) - {Dhu){tj) - hB,{tj/h + l/2){Fhu){tj)\ < \\Ehu\\Lo^ < C\\Ehu\\i < Ch\\u\\2. 
j 

However, B^{tj/h + 1/2) = 5i(l/2) = and therefore 

\\IhU - Dhu\\L^ < Ch\\u\\2 and \\Dhu\\L^ < C\\u\\2. (77) 



Theorem 3. Let iph be the solution of (29) with right-hand side g = Wip and e ^ 
Then 

max\Lpj - Lp{tj)\ = \\Lph- h'^h^ < C'/illv'lls- 

Proof. We rely on the first order asymptotic formula of Corollary [Tj Let Ch be the solution 
operator associated to (29), namely Ch^p = ^Ph- Let ^ := C_i{s)W~^{aip" — a^p). Then 
(72) shows that 

which can also be written as 

By the inf-sup condition in Theorem [TJ it follows that 

\\ChV-DHip-hChao<CeMs, 



-1 



and therefore, by (75) applied to C,, 

- DhV - hDh^o < h\\Ch^-Dh^\\o + Ch^Ms<Ch\U\\2 + \\iph) 

< C'h'M,, (78) 

since ||^||2 < C||9?||3 Note that for piecewise constant functions on a uniform grid of 
meshsize h we can estimate ||p/i||l°° ^ ^""'^''^IIP/illo- Thus, 

WCh^p - IhV\\L°° < h'^^/^WCh^p - Dh^p - hDh^\\o+ \\Dhip - Ih'^Wl^- + h\\Dh^\\L°° 

< ch^^^Ms + chyy + chUh, 



by (77). This proves the result. □ 



Theorem 4. Let (ph be the solution of (29) with right-hand side g = Wip and e G 
{-1/6,1/6}. Then 

max \'pj - 'p{tj)\ = \\(ph - 4v5||l°° < C/i^||v9||4. 

j 

Proof. The proof of this estimate is very similar to that of Theorem |3} We need to rely on 
the second order asymptotics of the error (Proposition [l2| to reveal the first non- vanishing 
term in the asymptotic error expansion when e G { — 1/6, 1/6}. In addition to this, using 
Proposition 14 and the Sobolev imbedding theorem, it is easy to show that the estimate 
(77) can be improved to \\IhU — DhuWioo < Ch'^\\u\\3. An inverse inequality, the stability 
estimate (Theorem [T]) and Theorem |2j can be used to show that 

\\(ph - Dhif - K^Dh'^Wioo < C/i^/^||v9||4, where 7 := W^V^^p, 



and P| G S{2) is given in Proposition 12 All remaining details are omitted. □ 
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7 Numerical experiments 



We will now illustrate some of the previous convergence estimates with a simple example. 
We take two ellipses, one centered at (0, 0) with semiaxes 1 and 2, and a second one 
centered at (4, 5) with semiaxes 2 and 1. We look for solutions of (18) (radiating solutions 
of the Helmholtz equations) in the exterior domain that lies outside both ellipses, with 
Neumann conditions on the boundaries (see Section 2.3). As exact solution we take 
U{z) := Hl)^\k\z - zol), where zq := (0.1, 0.2) is a point inside the first of the obstacles. 
We have taken A; = 1 in all examples. 



Experiment 7^1 (indirect method). After parametrization of the ellipses, a double 
layer potential (17) is defined on each of the curves. They are then used to set up a 2 x 2 
system of parametrized boundary integral equations, with diagonal terms of the form ([s]) 
and integral operators with smooth kernels as off-diagonal terms. We solve the system 
and plug the resulting densities in the fully discrete double layer potentials (20). We 
compute the errors: 



en := max |f/(z) - Uh{z)\, where O := {(4, 4), (5, 5.5), (6, 7), (7, 7.6), (6.8, 6)}, 
zeo 

that is, we observe the difference of the exact and discrete solutions at five external points. 
We expect Ch = 0{h) (this follows from Theoremj2]) and = 0{h'^) when e E {—1/6, 1/6} 
(Corollary [2]). The results are shown in Table [ij To see how the superconvergent values 
of e are reflected in the error, we plot the error Ch as a function of e for a fixed value of 
N in Figure [TJ 



N 


error 




e.c.r 


N 


error 




e.c.r 


10 


4.3005E(- 


-002) 




10 


9.7262E(- 


-003) 




20 


1.9193E(- 


-002) 


1.1640 


20 


2.5602^(- 


-003) 


1.8995 


40 


9.0917E(- 


-003) 


1.0779 


40 


6.2157E(- 


-004) 


2.0645 


80 


4.4279E(- 


-003) 


1.0379 


80 


1.5443E(- 


-004) 


2.0090 


160 


2.1852E(- 


-003) 


1.0189 


160 


3.8588^(- 


-005) 


2.0007 


320 


1.0855E(- 


-003) 


1.0094 


320 


9.6507^(- 


-006) 


1.9995 


640 


5.4097E(- 


-004) 


1.0047 


640 


2.4135^(- 


-006) 


1.9995 



Table 1: Errors and estimated convergence rates for Experiment #1. The variable is 
the number of points on each of the curves. The leftmost table corresponds to e = 1/3 
(order one) and the rightmost table to e = 1/6. 



Experiment ^2 (Richardson extrapolation). With the same geometric configura- 
tion, exact solution, and numerical scheme as in the superconvergent case (e = 1/6), we 
apply Richardson extrapolation to propose the potential 

Uh '■= tUh/2 - ^Uh, 
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0.2 0.4 0.6 0.8 



Figure 1: Error as a function of e in Experiment # 1. The superconvergent methods 
can be clearly seen as kinks in the error graph (corresponding to the first term in the 
asymptotic expansion of the error going through a zero). The methods becomes unstable 
as e — i- Z. Although our analysis does not cover this case, it is clear from the graph that 
e = 1/2 continues smoothly the error graph. 



as an improved approximation of the solution. The result of Proposition [12] points clearly 
to the existence of an asymptotic expansion of the error, very much in the style of those 
obtained for operator equations of zero or negative order in [2]. The numerical result 



shown in Table |2j corresponding maximum errors '■= maXzgo — = 0{h ). 

Experiment ^3 (direct method). We now apply a direct boundary integral equation 
method for the same exterior Neumann problem as in the previous experiments. This leads 
to a 2 X 2 system with the same matrix of operators as in the previous formulation, but the 
adjoint double layer operator appears in the right-hand side of the system. This operator 
is simply discretized with midpoint formulas on each of the intervals: see [1] for a similar 
treatment in systems related to the single layer potential. With this formulation, the 
unknown is the parametrized form of the trace of the exact solution = f/ o x and we 
can thus compare errors (Theorems [s] and |4]). We measure maximum absolute value 
of errors for ip on the points tj. The results are reported in Table [3] 

Experiment 7^4 (condition numbers). In this final experiment, we pick the matrix 
of the previous examples and compute its spectral condition number. We then show how 
a Calderon preconditioner based on premultiplying the matrix Wjj by a matrix |2], H] 

V., = ff«(/^|x(t.)-x(t,+,)|) (79) 

reduces the condition number of the resulting system to what is basically a constant 
/i— independent condition number. 
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N 


error 




e.c.r 


1 n 




UUU J 




20 


7.4235E(- 


-008) 


5.8707 


40 


5.6231E(- 


-009) 


3.7227 


80 


6.6107E(- 


-010) 


3.0885 


160 


8.1033E(- 


-Oil) 


3.0282 


320 


1.0052E(- 


-Oil) 


3.0110 



Table 2: Errors and estimated convergence rates for Experiment ^2: Richardson extrap- 
olation applied to one of the surperconvergent methods [e = 1/6). Errors are computed 
in several external observation points. The result reported with = 20 uses a grid of 20 
points as /i— grid and a refined grid of 40 points as /i/2— grid. 



N 


boundary 


error 


e.c.r 


N 


boundary 


error 


e.c.r 


10 


4.8555E(- 


-001) 




10 


2.0406E(- 


-001) 




20 


1.3426E(- 


-001) 


1.8546 


20 


3.5629E(- 


-002) 


2.5179 


40 


5.4891E(- 


-002) 


1.2904 


40 


8.6392E(- 


-003) 


2.0441 


80 


2.4641E(- 


-002) 


1.1555 


80 


2.1497^(- 


-003) 


2.0068 


160 


1.1792E(- 


-002) 


1.0632 


160 


5.3603^(- 


-004) 


2.0037 


320 


5.7677E(- 


-003) 


1.0318 


320 


1.3385E(- 


-004) 


2.0017 


640 


2.8527E(- 


-003) 


1.0157 


640 


3.3444E(- 


-005) 


2.0008 



Table 3: Errors and estimated convergence rates for Experiment ^^3. The leftmost table 
corresponds to e = 1/3 (order one) and the rightmost table to e = 1/6. The table shows 
errors \\ifh- hfHi^- 



A Second order asymptotics 



This section contains the proof of Proposition [12} Note that this result is required for the 
proof of L°° convergence of the super convergent methods. In order to prove Proposition 
12 we have to go one term further in the different asymptotic expansions that were used 



in the proofs of Propositions M and 10 



A.l Technical background 

Lemma 11. There exists C such that for all e and h 

Proof. It is based on the same ideas as the proof of Lemma [9} using the inequality 

/l3 



^ 2 



u{t)dt — hu{c) 
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«"(c) 



^ 2 



as starting point. 



□ 
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N 


V V V 


pon n W 

^_/VJlH-t V V 


lU 


/? nv^ AO 

0.9546 


r 'vo -1 o 


zu 




1 1 7QQ9 


40 


6.5349 


23.7403 


80 


6.5196 


47.5489 


160 


6.5159 


95.1320 


320 


6.5150 


190.2811 


640 


6.5148 


380.5709 



Table 4: Condition numbers for the matrix W of Experiments # 1 and # 2 and for the 



matrix VW, with V given by (79). 



Proposition 13. Let 01,02,03 G V and consider an operator 

V := aiHD_i + asHD.s + asHD^s + K4 where K4 G ^(-4). 

Let then Li := ai, L2 := o-iD — a2, L3 := aiD^ — 2a2D + 03, and consider the operators 

2 

R*u := Vm - Vg^^F^M + 5^ /i^C:,( ■ //i)L,F;,M, 

£=1 

3 

T*u := V«-VQ-l/2^^n + 5^/i^C:,(-//i + l/2)L,F;,«. 



Then 



fuWo + h\\R*u\\i < Ch^W'uh VueH'^, 



||rfM||o + /i||TfM||i < C/i^||m||3 ^u&h-\ 
Proof. It is a direct consequence of Proposition 16]. 

Proposition 14. Let E*u := u - DhU + Y.]=i ■ /h + l/2)FhU^^\ Then 

\\E*u\\q + h\\E*u\\i < Ch^\\u\\3 Vm G H\ 
Proof. See [5l Proposition 1] and the proof of [2, Proposition 16]. 



(80) 
(81) 

□ 



□ 



A. 2 Proof of Proposition 12 



Following (55) and (62), we consider the decomposition of the consistency error in five 
terms 

Wh{^h,Vh-Dh^p) = Wh{^h,'fh-^)+Wh{'iph,V- DhV) = (Ti + T2 + T3) + (T4 + T5). (82) 
To bound Ti we use Lemma [TT] with u = DViD: 

Ti - g{Q;:>/^,D3ViD<^}| < Ch^AWoMU- (83) 
25 



Proceeding as in (57) we can bound 

\T2\<Ch'UH\ 



own- 



To expand T3 we use Proposition 13 applied to V2 and (50) to obtain 



Tic 



We now use Lemma [9[b) and (50) to bound |T3a| + |T3;,| < C/i^H-i/^/iHoHv^lh, as well as 
Lemma gb) and (ISOl) to bound \T^c\ < Ch^Wi'hUWvh- Therefore 



IT3 + hC,{e){Q-',l^u, a^) + h''C^{e){Ql\xljH, (aD - b)^)\ < Ch^U^hM^. (85) 



To expand T4 we use Proposition 13 applied to Vi. Note that Li = al, L2 = aD — c and 



L3 G S{2). Because of (63), we can write 
2 



1=1 

2 



rnt 



h'C,{e){ijlUF,y,') + {ijlT*^'). 



— 

T4c 



Using Lemmata) and (48) we can bound 



\Tl\ + \Tl\ < Ch' 



'jhWo (llDLiDy^lh + llDLaDv^lli) < C'h^MolWU- 



By (50) (and using the commutation DFh = F^D to simplify some expressions) we next 
bound 



\Tl,\ + \Tl\ < ||7/>,||oU||DLiD(<^-F,vp)||o + /i2||DL2D(99-F;,<^)||, 



< cwijMhy-Fhiph + ti'y-FhM <c'h'\\iJh\\o 



Similarly 



Finally, by (81) 



\T,c\ < /l'|C:3(£)|||^;,||o||DL3DF;,(^||o < C || ||o || ^||4. 



T4d\ < ||^/.||o||DT#^'||o < CeU4o\MA. 



26 



Collecting all these bounds we have just proved that 



T4 -Y,h'C,{e){Q-^l^ilJh,'DU^^) 



< Ch^\\^hh\W\W 



^6) 



We are only left to deal wfth T5. Using Proposition 14, the argument in (69), and the 
fact that 52(1/2) = —1/24, we can write 

Q-^\^-Dr.^) = -hB,il/2)Q-,'F,if'-h'B,imQH'FHV" + Qk'E*^ 



24 



Therefore, 



h 



T5 = :^{QH!ei^k,y2Q-,'F,^"} + {Q,li;,,V2Q-,'E*^} 



24 
24 



where we have applied Proposition [7| By Lemma |9|^b) and (53) we can bound |T5a| + 



T5b\ < Ch' 







ly^lls, while by Lemma 10 and Proposition 14 we can bound \T^c\ < 



llV^felloll'y'lla- Therefore 



< Ch^\ 



11) 



The result is the combination of (82)-(87). 
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